Temperature-dependent "phason" elasticity in a random tiling quasicrystal 
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Both "phason" elastic constants have been measured from Monte Carlo simulations of a random- 
tiling icosahedral quasicrystal model with a Hamiltonian. The low-temperature limit approximates 
the "canonical-cell" tiling used to describe several real quasicrystals. The elastic constant K2 changes 
sign from positive to negative with decreasing temperature; in the "canonical-cell" limit, K2/K1 
appears to approach —0.7, about the critical value for a phason-mode modulation instability. We 
compare to the experiments on i-AlPdMn and i-AlCuFe. 

PACS numbers: 61.44.Br, 62.20.Dc, 5.10.Ln, 64.60.Cn 
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Quasicrystals possess translational long-range order, as 
manifested by resolution-limited Bragg peaksEI, yet pos- 
sess symmetries (e.g. icosahedral) incompatible with pe- 
riodicity. The best-understood quasicrystals can be rep- 
resented as rigid tilings (or ch^ter networks) decorated 
by atoms in a uniform fashiorojlffl. Their long-range or- 
der can then be investigated in a purely tilinfriramework, 
being determined by the "tile Hamiltonian" aa), which as- 
signs to every tiling the corresponding structural energy. 

Despite years of study, there remain two plausible, 
competing scenarios that predict quasicrystal long-range 
order (in the sense that both predict Bragg peaks in di- 
mension d > 2.) The "ideal tiling" scenario postulates 
that the tile Hamiltonian enforces an essentially unique, 
perfectly quasiperiodic ground state analogous to a p£ 
feet crystal (just as matching rules or covering rules 
enforce a Penrose tiling). Alternatively, the (equilibrium) 
"random tiling" scenario supposes an ensemble of many 
nearly degenerate packings of the structural units. Then 
the true ground state is normally a coexistence of crys- 
tal phases, but at higher temperatures the random-tiling 
quasicrystal becomes thermodynamically stable owing to 
its larger entropy (See references cited in |J) . 

The two scenarios may be distinguished experimentally 
by their diffuse scattering. In the random-tiling case, ev- 
ery Bragg peak is surrounded by wings of diffuse scat- 
tering with distinctive shapes predicted by the "phason" 
elastic theoryQ (See (||), below.) Those shapes were£>b» 
served in X-ray, and neutron diffraction of i-AlPdMn0'E2l 
and i-AlCuFeO, and the fitted phason-elastic constants 
were temperature-dependent. 

In this Letter, we simulate random tilings of rhombo- 
hedra with a tile Hamiltonian that favoxa_configurations 
approximating a "canonical cell tiling"BE3 (CCT). Our 
motives are (1) to infer the phason elastic constants of 
the CCT , which is difficult to simulate directly, and (2) 
to observe the crossover with decreasing temperature T 
as clusters of small (rhombohedral) tiles get bound to- 
gether into larger (CCT) "supertiles" in a toy model. 
The same techniques will be essential in computing the 
phason elastic constants for models of specific real qua- 
sicrystals. 



CANONICAL-CELL TILING AND 
HAMILTONIAN 



Atomic structure models of i- AlZnMgp [or i-AlM nSitj] , 
and many related alloys, may be built by placing icosa- 
hedral "Bergman" [or "MI"] clusters on the nodes of net- 
works with inter-node linkages of two lengths: b = 2.75a/? 
along twofold symmetry directions or c = 2.38a/? along 
threefold directions, which we shall call a "frc-network." 
(The quasilattice constant clr is set to unity in this pa- 
per; a# Ri,0.5nm in real quasicrystals.) A "canonical- 
cell" tilinga (CCT) is a special 6c-network built from the 
four smallest polyhedral cells having b or c edges. 

A Monte Carlo simulation of the CCT is practically 
intractable, since local tile reshufflings do not exist: spa- 
tially extended clusters must be rearranged to reach an- 
other valid tilingBo. To circumvent this problem, our 
model system is a be network in which the CCT con- 
straints are not imposed, but are favored by a Hamilto- 
nian and thus satisfied in the limit T — > 0.— We emphasize 
that unlike matching-rule Hamiltonians,o eq. ([!]) does 
not force a unique quasiperiodic T = state; the ground 
state ensemble is still a random tiling, but a different one 
(with less entropy) than at T = 00. Our configuration 
space is the random "three-dimensional Penrose tiling" 
(3DPT), which includes all packings of the well-known 
oblate and prolate rhombohedra having edges of length 
an = 1 along 5-fold symmetry directions. Similar simu- 
lations of another imperfect approximation of the CCT 
by a bc-network, namely the lattice-gas of Ref. [l2| will 
be reported elsewhere. EJ3 

The 3DPT contains a special subset of vertices, called 
"12-fold" nodes, at which rhombohedron edges emanate 
in all 12 of the five-fold symmetry directionsEB. This 12- 
fold network approximates the CCT, since the separation 
of nearby nodes is usually a b or c linkage; however a 
"short" separation of or (along 5- fold axes) occurs rarely. 
For our simulations, we adopted the Hamiltonian 



Hu = -Nl2 + ashort-ZVsb, 



(1) 



where N12 is the total number of 12-fold nodes, and N s b 
is the number of "short" node pairs with separation or. 

The purpose of cohort is to suppress "short" pairs; 
cohort = 1 suffices, which we adopted after testing other 
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values. We conjecture (and simulations support) that the 
ground states of (Q) are a subset of CCT's, namely those 
with maximum node density, and not the maximum- 
entropy random CCT. 



II. PERPENDICULAR SPACE AND 
FLUCTUATIONS 



Our rhombohedron edges are ±e'l a , whereO 



(r, 0, 1)/(1 + t 2 ) 1 / 2 and the other basis vectors are made 
by reflections in the x, y, or z planes, or by cyclic per- 
mutation of the components. Thus every vertex has co- 
ordinates of the form 



(2) 



In a standard trick, the integer coefficients may be visu- 
alized as coordinates x = [m, . . . , uq] of a 6-dimensional 
hypercubic lattice, into which eq. (||) embeds a 3- 
surfacala. This 3-surface is conveniently parametrized by 
a "perpendicular" ("perp") coordinate, defined for each 
vertex by 



a=l 



e-L 



(3) 



where e ± i = ry(l,0, — r)/(l + r 2 ) 1 / 2 and the other basis 
vectors e a are derived by the same reflections and cyclic 
permutations used to produce e" Q from e"i. We write 
h(r) for the smoothed (coarse-grained) version of x^. 

The gradients Vh(r), forming the "phason strain" ten- 
sor, quantify the local deviation fcpm icosahedral sym- 
metry. The random-tiling scenarioel predicts that a qua=. 
sicrystal's free energy density has a phason elastic form,t£l 
i.e. is quadratic in components of Vh. This was con- 
firmed numcijically_fjor several completely random (T = 
oo) tilingsBBElM'El, but not until now for the CCT. In 
contrast, the "ideal tiling" behavior implies a free en- 
ergy scaling as |Vh|. An icosahedral quasicrystal has 
two elastic constants, K\ and K%. Very roughly speak- 
ing, K\ parametrizes the mean strength of the elastic 
stiffness, and parametrizes its anisotropy (degree of 
coupling between direction of the gradient in real space 
and components of h). 

For fluctuations around a state with zero mean phason 
strain, the dimensionless elastic free energy isEZlEj 



F / T = j£C y -(q)Mq)M-q) 



(4) 



wehere it is convenient to use h(q) 

y-l/2 Jrf3 re -iqT h ( r ) Here 

Cy(q) = KilctfSa-Ki [(i|q| 2 +2 (Z 2 
+Qi+jT - rql_ x ) 8 t] - 2q iqj ] . 



Eq. (H) implies that each h(q) consists of Gaussian ran- 
dom variables, and their correlations 



T y (q) = (^(-q)^(q)). 



are given by 



r«(q) = [C(q)- 1 ] 



(0) 



(7) 



The diffuse diffraction intensity^ near a Bragg peak with 
reciprocal lattice vector G" is 



(8) 



where G 1 - is the perp-space partner of G". (An ecu 
librium random tiling in d > 2 has true Bragg peaksa 
as well as diffuse scattering.) Experiments measure 
simulations can measure <M) directly. 



III. SIMULATION, MEASUREMENTS, AND 
CORRECTIONS 

Our data is taken from tilings using periodic boundary 
conditions in a cubic cell of side L = t 6 (2— 2/ v^) 1 ^ 2 ^ = 
18.868afl. This corresponds to the "8/5 Fibonacci ap- 
proximant" of the 3DPlfe3, and contains N v = 10336 
rhombohedron vertices. We use the very simple update 
move for the 3DPT which rearranges the four tiles sur- 
rounding a vertex with just four edges .ljlj. 

Simulation time is measured in Monte Carlo sweeps, 
where each sweep makes N v actual flips on randomly cho- 
sen vertices. A simulation run consisted of 500 thermal 
cycles; each cycle included 3000 MCS at T = oo (to erase 
memory of the preceding cycle), followed by 34 cooling 
stages each at a constant temperature. Each stage takes 
1000 MCS, except the first six stages (T > 0.8) used 200 
to 800 MCS; 10 measurements were made during the last 
100 MCS of each stage. A typical 8/5 cooling run took 
10 days of CPU time on a 1.5 GHz Athlon processor. 

We constructed x" and x 1 - for each vertex, using (j|) 
and (|J). For each q on a standard list of 70 wavevectors, 
we took 



v x£3DPT 



x -L e -iq-x" 



(9) 



(5) 



as in Rcf. g2| where V is the simulation cell volume 
and the sum runs over the N v actual 3DPT vertices. 
The normalization in (^) ensures that H(q) ~ h(q) 
for small q. Fluctuation correlations analogous to (|g), 
r H jj(q) = (Hj(— q)Hj(q)) were accumulated, and were 
symmetrized with respect to cyclic permutations (xyz). 

Although (h(q)) = 0, Eq. (§) has a nonzero expec- 
tation when q is a Bragg vector, just like the structure 
factor, from which the sum (^|) differs only by the extra 
factor x- 1 -. As in the structure factor, each Bragg peak of 
our r H (q) function is surrounded by a tail of diffuse scat- 
tering, which must also be subtracted in order to correct 
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the nearby wavevectors. To calculate this tail, we adopt 
the "undulating cut approximation" (Rcf. ||, Sec. 7), 
in which the occupation probability of a physical-space 
vertex x" is /9o( x_L — h( x )) < 1- Here x" and x x arc 
related as in eqs. (^) and (||), h(r) is smooth and fluctu- 
ating according to (0), and po(h) defines an ideal refer- 
ence structure when h is flat. The measured fluctuation 
is predicted to be&J 



r H ii(q) 



r y (q) + E r r gs (q> G ' 



(10) 



G 



where 



r Bra gg(q;G || )=A ^ 



vs. 



q,G J 



(G x ,r(q-GU)G J 



Here D = D(G _L ) can be expressed^ in terms of the 
Fourier transform of po(x ); in light of the symmetry, 
D(G ) || G along any symmetry axis. 

Our measurements include one orbit of weak Bragg 
peaks on 3-fold symmetry axes, - at G" = ^(333) and 
^(520) in our cell - where the Bragg intensity con- 
tributes far stronger than the long-wavelength fluctua- 
tions (H). We determined each Bragg term (first term in 
(|TT|) , i.e. AA™) from (h(q)) averaged over each tem- 
perature stagecJ. ^From our estimate of K\ and K 2 , 
we know T(q) via (Q), hence the second term in (11), 
so we can subtract the sum in (|l0|) from the measured 
data. The resulting matrices T(q) are inverted to obtain 
(via (0)) estimates of (||), which are fitted to K\ and K 2 
using linear least squares. The process is iterated until 
converged. 



IV. RESULTS AND DISCUSSION 

To clarify two other artifacts which depend on |q|, 
we fitted effective constants il Ki(n)" and u K 2 (n) v sep- 
arately for each shell n defined by |q| 2 = n(27r/L) 2 ; the 
results are shown in Fig. |], for selected temperatures. 
The continuum approximation must break down at the 
CCT linkage length scale (~ 2.5<zr), which explains the 
deviations visible forn > 20, i.e. 27r/|q| < 4.2a#. (Above 
the supertile formation temperature, the tile length scale 
is the rhombohedron edge o_r, and no such deviations are 
seen in this wavevector range.) 

At small wavevectors, equilibration time becomes 
problematic. As shown in Ref. eqs. (20ff), the re- 
laxation rates of h(q) are 76^ (q), where 7(T) is a re- 
laxation coefficient and {Ci, 62,63} are the eigenvalues 
of C(q). This becomes problematic when temperature 
is low (since 'j(T) — ► as T — > 0) and |q| is small (as 
C cx |q| 2 , eq. (g). The consequence is that, at the lowest 
temperatures in Fig. |l|, the smallest wavevectors are gov- 
erned by the elastic constants of the higher temperature 
at which these modes fell out of equilibrium. 

The plots in Fig. @ suggest two-regimes as a function 
of wavevector and temperature.Eil At high T or large 



wavevector, a linear dependence on n oc |q| 2 is visible, 
with the slope increasing as T — ► 0. We attribute this 
to our use of the 3DPT vertices in (||): the x 1 - differ- 
ences between adjacent vertices (separated by an) cer- 
tainly exceed the fluctuations of the smooth h(r) sur- 
face, extrapolated to this short scale. Consequently, the 
inferred elastic constant is spuriously small. 

On the other hand, a second regime is visible at low 
T and |q| < k(T), where k(T) defines a crossover. In 
this regime, Ki(n) has a clear linear dependence, which 
we tentatively attribute to higher-order-gradient terms 
of 0{q A ) in the harmonic phason-elastic coefficients (||), 
which are permitted by symmetry. Although it certainly 
affects the n = 1, 2 wavevector shells at the lowest tem- 
peratures, lack of equilibration cannot explain all of the 
low T/small |q| regime, as the data (in this.regime) agree 
very well with our lattice-gas simulationsE3 which have a 
different dynamics. We used the empirical fitting form 



Ki(n) 



\Jn 2 + ■ 



(12) 



to extrapolate Ki(n) and Ki(n) to q = 0. 

Our principal result is the temperature dependence of 
the elastic constants indicated on Fig. ^. At infinite tem- 
perature, i.e. for the maximally random tiling, we fit 
Ki = 0.84 and K 2 = 0.52, so K2/K1 = 0.62 in agree- 
mept. with previous measurements on the same ensem- 
blcOEj. As temperature decreases, K\ increases, as ex- 
pected since the Hamiltonian (|l|) favors formation of "su- 
pertiles" which have small perp-space differences among 
their corner nodes. Meanwhile, K2/K1 turns negative 
with decreasing temperature_-ponsistent with our knowl- 
edge that the random CCTEa has K 2 < 0. It appears 
that Ki 2.26 and K 2 -> -1.44 as T -> 0, thus 
K 2 /K\ — ► —0.64- with errors of order 10%. Our lattice- 
gas simulationalj gave similar results. 

Approaching T = 0, the fit of our data to the angular 
dependence of elastic theory may be worsening. Further- 
more, K 2 jK\ seems to be approaching the critical value 
—0.75, at which the system goes unstable to a modula- 
tion with q and h(q) along a 3-fold axis. This suggests 
that the true canonical-cell tiling - more precisely the 
maximum-density subensemble of the CCT - might have 
a more complicated dependence than elastic theory. This 
is not unprecedented: exactly that is known to happen 
in those. 10-fold and 12-fold symmetric two-dimensional 
tilingstj that have extra constraints like the CCT. 

Neutron and X-ray diffraction measurements of diffuse 
wings amuni. Bragg peaks yield K 2 jK\ — —0.5 for i- 
AlPdMrMi^y, agreeing with our result, so i-AlPdMn 
might be modeled by a CCT-like network of microscopic 
clusters. On the other hand, i-AlCuFe has K 2 jK\ > 
fjlil. Thus, despite their similar face-centered icosahe- 
dral atomic structures, z-AlPdMn and i-AlCuFe must be 
modeled by different tile Hamiltonians. 

The recent experimental calibration of i-AlPdMn 
dataS yielded K x /T = 0.002A~ 3 in absolute physical 
units. Multiplying by Oij 3 (with an — 4.5A) converts this 
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to K\ w 0.2 in the units of this paper, which is puzzlingly 
small. This is hard to explain by a misjudgment in our 
model formulation. If our tiles are invalid, one expects 
that (since this alloy is so well ordered), they should be 
replaced by large supertiles, so that K\ becomes larger, 
or that matching rules are satisfied implying K\ — ► oo. 
Instead, the data suggest that, implausibly, atoms are 
free to fluctuate in groups even smaller than a rhombo- 
hedron. 

The temperature dependence of (T) was used for ex- 
trapolation only; the Hamiltonian (|lj) simply implements 
the CCT and is unlikely to model a realistic quasicrys- 
tal. A physically relevant model would add to ([!]) an 
additional perturbation a' that breaks the degeneracy 
among canonical-cell tilings, so that large a' corresponds 
to physical low temperature. Such a model can be inves- 
tigated using the techniques of this paper. 



In conclusion, we evaluated the phason elastic con- 
stants of the canonical-cell tiling for the first time, by a 
simulation in which the mean-square fluctuations of the 
abstract surface representing the tile configuration were 
measured at many wavevectors. This is the first the- 
oretical measurement in any quasicrystal model of the 
temperature dependence of the phason elastic constants. 
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FIG. 1: Phason elastic constants K\ (bottom) and K2/K1 
(top), where K 1 and K2 are fitted from "perp" space fluctua- 
tion r H ; m (q) of rhombohedron vertices, Fourier-transformed 
and separately fitted in each shell of wavevector q = 
■Jntliir/L). Bragg peaks and tails were subtracted using 
(p.0|). Only five of the 34 temperatures are shown (labeled 
by P = 1/T. The fits to Eq. (113) are shown as solid lines. 




FIG. 2: Elastic constants as a function of inverse tempera- 
ture 8. These are same fitted values indicated by solid lines in 
Fig. qL extracted from the perp-coordinates of rhombohedron 
vertices and extrapolated to q — > 0. 



